# Title     : TODO
# Objective : TODO
# Created by: Administrator
# Created on: 2019/7/24

library(ggrepel)
library(optparse)
library(magrittr)
library(ggpubr)
library(tidyverse)

createWhenNoExist <- function(f) {
  !dir.exists(f) && dir.create(f)
}

option_list <- list(
  make_option("--i", default = "AllMet.csv", type = "character", help = "metabolite data file"),
  make_option("--g", default = "SampleInfo.csv", type = "character", help = "sample group file"),
  make_option("--sc", default = "sample_color.txt", type = "character", help = "sample color file"),
  make_option("--config", default = "config.csv", type = "character", help = "config file")
)
opt <- parse_args(OptionParser(option_list = option_list))

configData <- read.csv(opt$config, header = F, stringsAsFactors = F) %>%
  set_colnames(c("arg", "value")) %>%
  column_to_rownames("arg")

configData

diffMethod <- configData["diffMethod.mulMethod", "value"]

my.test <- function(data, method) {
  p <- 0
  curMethod <- method
  if (curMethod == "auto") {
    sTest <- data %>%
      group_by(ClassNote) %>%
      summarise(p = {
        var <- var(value)
        if (var == 0) {
          0
        }else shapiro.test(value)$p.value
      })
    bTest <- bartlett.test(value ~ ClassNote, data)
    bp <- bTest$p.value
    if (all(sTest$p > 0.05) && bp > 0.05) {
      curMethod <- "anova"
    }else {
      curMethod <- "kw"
    }
  }

  if (curMethod == "kw") {
    test <- kruskal.test(value ~ ClassNote, data = data)
    p <- test$p.value
    if (is.na(p)) {
      p <- 1
    }
    curMethod <- "kruskal.test"
  }else if (curMethod == "anova") {
    test <- compare_means(value ~ ClassNote, data = data, method = "anova")
    p <- test$p
    curMethod <- "anova"
  }
  list(p = p, method = curMethod)
}


options(digits = 3)

sampleInfo <- read.csv(opt$g, header = T, stringsAsFactors = F) %>%
  select(c("SampleID", "ClassNote")) %>%
  mutate(ClassNote = factor(ClassNote, levels = unique(ClassNote)))
sampleIds <- sampleInfo$SampleID

groups <- unique(sampleInfo$ClassNote)

allData <- read.csv(opt$i, header = T, stringsAsFactors = F) %>%
  select(-c("HMDB", "KEGG")) %>%
  select(c("Class", "Metabolite", sampleIds)) %>%
  rowwise() %>%
  do({
       result <- as.data.frame(.)
       kwData <- result %>%
         select(-c("Metabolite", "Class")) %>%
         gather("SampleID", "value") %>%
         inner_join(sampleInfo, by = c("SampleID"))
       for (group in groups) {
         sampleIds <- sampleInfo %>%
           filter(ClassNote == group) %>%
           .$SampleID
         groupData <- result[sampleIds] %>%
           unlist()
         result[, paste0(group, ".Mean")] <- mean(groupData)
         result[, paste0(group, ".Median")] <- median(groupData)
         result[, paste0(group, ".SD")] <- sd(groupData)
         iqr11 <- quantile(groupData, 0.25)
         iqr12 <- quantile(groupData, 0.75)
         iqr1 <- paste0("[", iqr11, ",", iqr12, "]")
         result[, paste0(group, ".IQR")] <- iqr1
       }

       rs <- my.test(kwData, diffMethod)
       result$P <- rs$p
       result$test.method <- rs$method
       result
     }) %>%
  select(-c(sampleInfo$SampleID)) %>%
  as.data.frame()

data <- allData %>%
  mutate(FDR = p.adjust(P, method = "fdr")) %>%
  select(-"test.method", everything(), "test.method")

print(head(data))

parent <- paste0("./")
allFileName <- paste0(parent, "/AllMet_Test.csv")
write.csv(data, allFileName, row.names = FALSE)












